goodness gracious
big balls of HAIL

Use R, leaflet.js, and dygraphs.js to assess hail risk

Author 1, Affiliation
Author 2, Affiliation



Imagine baseball-sized ice balls falling from the sky.

For some Americans, this is a real weather phenomenon that risks the physical well-being of people and property. But how often does it happen? Where do these events normally happen?

As it turns out, the National Oceanographic and Atmospheric Administration (NOAA) collects massive amounts of data on precipitation using radar stations across the country. Algorithms can be built on top of the data to track ‘hail signatures’ detecting when these events occur and the severity. This severe storm data is captured in NOAA's Severe Weather Data Inventory (SWDI) housed within the National Centers for Environmental Information (NCEI).

Looking at nearly 10 million events over the last 10 years, we can see that most of these events happen in the Midwest.

Hail Risk (2005 to 2015)

Daily Probability of Hail by 1/4 degree grids (~17 miles)

A closer look at severe events (e.g. hail diameter > 3in), the majority of these events occur in a small snaking section of the country indicating that it's relatively isolated.

Severe Hail Risk (2005 to 2015)

Daily Probability of Hail (Diameter > 3in) by 1/3 degree grids (~23 miles)



[Section to be completed] Get Started with SWDI
SWDI data can be used for everything from agriculture to risk management. 

setwd("/Users/jeffreychen/Google Drive/DOC/027-SWDI/cviz_tutorial")
    library(sqldf)
    library(RColorBrewer)
    library(leaflet)
    library(googleVis)
    library(ggplot2)
    library(DT)
    library(rgdal)
    library(sp)

    ##Set Parameters
    start_date = "2005-01-01"
    end_date = "2014-12-31"
    range <- as.Date(end_date,"%Y-%m-%d")-as.Date(start_date,"%Y-%m-%d")
    series = "nx3hail"
    fraction = 0.25

    #Load in files
    shp  <- readOGR("cb_2014_us_county_20m", "cb_2014_us_county_20m", stringsAsFactors=FALSE, encoding="UTF-8")
    raw <- read.csv("hail_10yr.csv")
    fips <- read.delim("http://www2.census.gov/geo/docs/reference/state.txt",sep="|")
    fips <- fips[,c("STATE","STUSAB")]
    fips$STATE <- as.numeric(fips$STATE)
    fips$STUSAB <- as.character(fips$STUSAB)


    #Cutdown
    raw <- raw[raw$LON<(-50) & raw$LON>(-140) & raw$LAT > 25,]

    raw$LON <- round(raw$LON/fraction)*fraction
    raw$LAT <- round(raw$LAT/fraction)*fraction

    ##De-duplicate by day, latitude and longitude
    deduped_day <- sqldf("SELECT DATE, LON, LAT, MAXSIZE
                         FROM raw
                         GROUP BY DATE, LON, LAT, MAXSIZE")

    deduped_day$lvl_2in<-0
    deduped_day$lvl_2in[deduped_day$MAXSIZE>2]<-1
    deduped_day$lvl_3in<-0
    deduped_day$lvl_3in[deduped_day$MAXSIZE>3]<-1


    ##Daily gridded frequencies
    singles <- sqldf("SELECT LON, LAT,COUNT(date) cnt, SUM(lvl_2in) cnt_2, SUM(lvl_3in) cnt_3
                     FROM deduped_day
                     GROUP BY LON, LAT")

    for(i in 3:ncol(singles)){
      singles[,i]<-singles[,i]/as.numeric(range)
    }
    write.csv(singles,"singles.csv",row.names=F)

    ##Set up spatial
    points<-singles
    coordinates(points)=~LON+LAT
    proj4string(points)=CRS("+proj=longlat +ellps=WGS84 +datum=WGS84")
    shp <- spTransform(shp, CRS("+proj=longlat +ellps=WGS84 +datum=WGS84") )

    ##Spatial join
    a<-over(points,shp)
    singles <- cbind(singles,a[,c("STATEFP","NAME")])
    singles$NAME <- as.character(singles$NAME)

    #Add state abbreviation via FIPS
    singles$STATEFP <-as.numeric(singles$STATEFP)
    singles <- merge(singles,fips,by.x="STATEFP",by.y="STATE",all.x=T, sort=F)
    singles$loc <- paste("Nearest County: ",singles$NAME,", ",singles$STUSAB,sep="")
    singles$loc[is.na(singles$NAME)]<-""
    singles <- singles[!is.na(singles$STUSAB),]

    #Subset data into separate frames (needed for layers in leaflet.js)
    in3 <- singles[singles$cnt_3>0,]

    #Vectorize data for inclusion for popup text
    #All points
    county_all <- singles$loc
    x_all <- singles$LON
    y_all <- singles$LAT
    pr <- round(100*singles$cnt,3)

    #Hail balls > 3in
    county_3 <- in3$loc
    x3 <- in3$LON
    y3 <- in3$LAT
    pr3 <- round(100*in3$cnt_3,3)

    #Popup
    content_all <- paste("<h3>",county_all,"</h3>Grid Point: ",x_all,", ",y_all,"<p>Prob of Any Hail <span style='color:red'><strong>: ",pr,"%</strong></span></p>")
    content_3 <- paste0("<h3>",county_3,"</h3>Grid Point: ",x3,", ",y3,"<p>Prob of Hail (> 3in)<span style='color:red'><strong>: ",pr,"%</strong></span></p>")

    pal <- colorNumeric(
      palette = "Set1",
      domain = pr
    )

    pal2 <- colorNumeric(
      palette = "Set1",
      domain = pr3
    )

    leaflet(width="100%") %>%
      setView(lat = mean(y_all), lng = mean(x_all),4) %>%
      addTiles('http://{s}.basemaps.cartocdn.com/dark_all/{z}/{x}/{y}.png') %>%
      addCircleMarkers(data = singles, lat = ~ LAT, lng = ~ LON,radius=(pr/5),
                       fillOpacity = 0.8,stroke = FALSE,
                       color = ~pal(pr), popup = content_all) %>%
      addLegend("bottomright", pal = pal, values = pr,
                title = "Prob(Hail)",labFormat = labelFormat(suffix = "%")
      )



    leaflet(width="100%") %>%
      setView(lat = mean(y_all), lng = mean(x_all),4) %>%
      addTiles('http://{s}.basemaps.cartocdn.com/dark_all/{z}/{x}/{y}.png') %>%
      addCircleMarkers(data = in3, lat = ~ LAT, lng = ~ LON, radius=2*(pr3),
                       fillOpacity = 0.8,stroke = FALSE,
                       color = ~pal2(pr3), popup = content_all) %>%
      addLegend("bottomright", pal = pal2, values = pr3,
                title = "Prob(Hail diameter > 3in)",labFormat = labelFormat(suffix = "%")
      )

Commerce Data Service

github